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We study the numerical solutions of the Dicke Hamiltonian, which describes a system of many 
two level atoms interacting with a monochromatic radiation field into a cavity. The Dicke model is 
an example of a quantum collective behavior which shows superradiant quantum phase transitions 
in the thermodynamic limit. Results obtained employing two different bases are compared. Both of 
them use the pseudospin basis to describe the atomic states. For the photon states we use in one case 
Fock states, while in the other case we use a basis built over a particular coherent state, associated 
to each atomic state. It is shown that, when the number of atoms increases, the description of the 
ground state of the system in the superradiant phase requires an equivalent number of photons to 
be included. This imposes a strong limit to the states that can be calculated using Fock states, 
while the dimensionality needed to obtain convergent results in the other basis decreases when the 
atomic number increases, allowing calculations that are very difficult in the Fock basis. Naturally, it 
reduces also the computing time, economizing computing resources. We show results for the energy, 
^ ■ the photon number and the number of excited atoms, for the ground and the first excited state. 



o 



< 



i 

i- 
C 



> 



I. INTRODUCTION 



The Dicke model describes a system of TV two-level atoms interacting with a single monochromatic electromagnetic 
radiation mode. One of this model's interesting properties is its quantum phase transition in the thermodynamic limit: 
when the coupling constant between the atoms and the radiations increases, the system undergoes a phase transition 
between the normal and the superradiant phases. This transition is an example of quantum collective behavior. In 
recent years the Dicke model has attracted renewed interest in the study of the phase transitions, entanglement in many 
q-i body systems, chaos in quantum mechanics, and applications in quantum optics and condensed matter physics. At 
the light of experimental realizations of the superradiant quantum phase transitions [l| , exploring numerical solutions 
to the model becomes of practical interest. The purpose of this work is to present numerical solutions which make 
easy solving the system, given that, for finite TV, only in the limits in which the interaction strength between the 
atoms and the electromagnetic field becomes null or infinit the system is exactly solvable. In particular, we employ 
coherent states [1| to improve the efficiency in obtaining numerical solutions of the Hamiltonian. Expectation values 
t ■ of observables of interest are also presented. 
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f"^ ' A. Description of atomic-field interaction 

t-H \ 

The interaction between TV two-level atoms and a multimode radiation field can be described in its more general 
form by a Hamiltonian formed by three terms: one associated to the radiation modes, one related to the atomic part 
and a third one which describes the interaction between themQ. 



03 ; H = Y / ^ka{a k + AJ z + hY / A=(4+ak)(J+ + J-) (1) 
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Each radiation mode has a frequency ujk and it is weighted by the number operator a\ak of this mode. For the atomic 
part A is the transition frequency, meanwhile J z , J+, J_ are pseudospin collective atomic operators between the 
ground state \b) and the excited state \a), which obey the SU(2) algebra; they describe the two-level atoms ensemble 
with an eigenvalue j of length TV/2. The J z = \a)(a\ — \b)(b\ operator gives the number of atoms in the excited state. 
The J + = \a)(b\ operator excites the atom and the J_ = \b)(a\ operator lowers it to the ground state. Finally, the 
interaction parameter 7^ depends principally on the atomic dipolar moment. The description in this work considers 
a single mode of radiation and it is called the Dicke model. Therefore, the Dicke Hamiltonian can be written as [5j : 

H Di cke = wa'fl + AJ Z + — L (a f + a) ( J+ + J_) (2) 
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There are some well-known approximations of this model. In the so called rotating-wave approximation the terms 
ak J- and a k J+, which represent interactions beyond the exchange of energy between non-perturbed atoms and the 
field, are neglected. In this approximation it is called the Tavis-Cummings model [a]. If we consider the simplest case 
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we have only one atom interacting with a single radiation mode. This is called the Jaynes-Cummings model Q- In 
this case the J z , J + , J_ operators are commonly written as the a zi a +1 cr_ operators, which follow the same algebra, 
but describes only a single two-level atom. 

Returning to the Dicke Hamiltonian, the interaction parameter is scaled by yN factor to have the proper thermo- 
dynamical limit. The Dicke model is exactly solvable in two limits, discussed in the next section. Analytical solutions 
which are exact in the thermodynamic limit have been found both for the normal Q and the superradiant phases. 
In general, for finite N there are no exact solutions. For this reason it is worth to explore numerical solutions. As it 
was mentioned above, the interest in solving the Dicke model for a finite iV comes from the fact that many systems 
could be described using this Hamiltonian and also it could help to understand more about the phase transitions in 
quantum systems. 



II. NUMERICAL SOLUTIONS TO THE FINITE DICKE MODEL 



A. The limit 7 — > 0: Solutions using Fock states 



At zero interactions, the Dicke Hamiltonian eigenstates are the tensor product between Fock states for the radiation 
modes and angular momentum eigenstates for the atomic part \^) = \n)\j, to), with j = N/2. This is possible because 
the atomic operators are pseudospin ones, as it was mentioned above. The energy spectra in this case is 

^7->0 = ojn + 2^ TO - (3) 

For finite values of 7 the interaction term mixes states with different number of photons n and different occupation 
numbers in the two level atoms j + to. To diagonalize the Hamiltonian by numerical means it is necessary to truncate 
the number of photons in the field. The Hamiltonian matrix elements are very simple: 



(n ; j, m \H\n; j, to) = (nu + mii)8 n > n 8 m > m 

+ Vn ^ n + 1 Sn '> n+1 + ^ 5n '> n ~ l ) ( 4 ) 

x (yj(J + 1) -to(to + 1) 8 m > tm+1 + + -m(m- 1) S m > rn _^ 

For each eigenstate there exists a number of photons which is enough to ensure the convergence of the solution. So, 
this value will be the one in which the energy will converge to a correct value and we call this number for the ground 
state n m ax- The dimension of the Hamiltonian matrix to be diagonalized is dimHFock = ((n m ax + l)(2j + 1)) 2 . As j is 
increased, the number of photons necessary to obtain convergence increases too, turning impracticable the numerical 
evaluation of the eigenstates. 

The Dicke Hamiltonian has a dynamical symmetry associated with the projectors of the symmetric and antisym- 
metric representations of the cyclic group C2, given by [2j, 

Ps = \(l + e™ k ), P A = 1 -{l-e™ k ). 

This symmetry allows the classification of the eigenstates in terms of the parity of the eigenvalues A = j + to + n of 
the excitation number operator 

A= \JP + 1/4-1/2 + ^ + ^0. 

Under the unitary transformation U = e ,7rA ; it is immediate that 

UaU^ = ae-™ = -a, U J+U^ = J+ = -J+, U J Z U^ = J z (5) 

with the corresponding hermitean conjugated expressions. 

It must be remarked that the eigenstates of the Hamiltonian must be simultaneous eigenstates of the parity operators 
Ps or P A . 
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B. The limit A — > 0: Solutions using coherent states 



After a rotation of the pseudospin operators by \ around the y axis, the J z operator becomes — J' x = — (J' + + J'_) 
and J x becomes J' z . This leaves the Hamiltonian as: 



H = U(Ja-~A(J' + + J'_) + -jL{(J + a)J' z . (6) 



2 + 

Then, we build a new annihilation operator, by shifting the original one. This gives us: 



A = a+^=J' z =a + GJ' z . (7) 



By substituting A in the Hamiltonian we obtain: 

H = u J (A^A-G 2 J / z 2 )-^(J' + + J'_). (8) 
This Hamiltonian has the same symmetries of the original one. 

Using the operators A, A^ a new basis is built, employing coherent states which are functions of the parameter G 
and the eigenvalues of the J z operator. We can obtain this states by applying powers of the shifted creation operator 
acting over the A vacuum state. 



OO yi ATI 

« = -Gm)= e -H 2 / 2 y^4-|0) h (9) 



Now, the states that describe the radiation are a mix of the field and the atomic states, and the interaction is described 
in terms of G = J'Jjj ms tead of 7. It would seem that the operators are not coupled, however, this is not true because 
the operator A contains both a and J z . The matrix's elements in this coherent basis are @: 

b (n';j,m'\H\n;j,m) b = u6 n < J m > . m {n - G 2 m 2 ) 

-A(Vj(i + 1) - m ( m + 1) s m ', m +ib( n ; m \n;m + l) b 
+ Vj(j + 1) - vnivn - 1) 5 m > m _ lb (n ;m \n;m- l) b ) (10) 
The cost we pay for using this basis is that it is not orthogonal anymore, so we have overlaps between states: 

min(n ,n) 



~ V" Vn'lnl f 2 , 

(n - k)\(n - k)W. 



Z_ - .. ,.,„.. MILI (-1) B " G ™+»-^, (11) 

fc=0 



G 

e 2 



min(n ,n) 



(n'-fc)!(n-Jfc)!fe! (_1) & ' ' 



k=0 



Notice that, if we denote by n' the eigenvalues of A^ A and by m! the eigenvalues of J' z , in limit A — > the 
Hamiltonian is exactly solvable in the new basis, with energies 



E A ^ = uj(n' -G 2 m' 2 ). (13) 

In this limit the ground state is doubly degenerate, with n' = and to' 2 = j 2 . The eigenstates with the appropriate 
C*2 symmetry arc 

|A ^ 0,±) = -±= {| - Gj)\jj) ± \Gj)\j j}} 

When A ^ for 7 >> uj and 7 >> A, the maximum value of n', which we will call N max , will be very small. 
The dimension of the Hamiltonian matrix to be diagonalized in this case is dimH^ = {{N max + l)(2j + l)) 2 , which 
can be orders of magnitude smaller than dimlip 0( ± , allowing for calculations that in the other basis would have been 
impossible. 
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Figure 1: Fock n max as function of 7 for values of j from 1 to 10. 
III. RESULTS 

Numerical solutions are presented for the ground and first excited states in resonance, A = u>, for both basis. We 
have calculated n max (N max ), and estimated the computing time in seconds, for several j in resonance employing the 
ground state; these results are shown in table 1. The energies are measured in units of u. The source code was written 
in Wolfram Mathematica. It allows to find the eigenvalues and eigenvectors of the Hamiltonian, and to evaluate the 
expectation values of all observables of interest. The computing time is calculated through the Mathematica command 
"Timing", with the code running in a computer with a processor AMD Turion 64 X2 at 2.00GHz. While it is clear 
that the absolute times reported here strongly depend on the platform employed to run the code, the trends and the 
relative values should be mainly independent of them. 

In Fig. 1 we show n max as a function of 7 for several j for the Fock basis, while in Fig. 2 we plot N max in the 
coherent state basis. It must be remarked here that the Fock n max is the maximum value that the a) a eigenvalues 
can take, while the coherent N max is the same for the A' A eigenvalues, so the Fock n max defines the truncation in 
the photon number and the coherent n max does not. It is worth to compare them because both are proportional to 
the dimensionality of the Hamiltonian matrices in each basis. As it is shown in the table and in the figures, for a 
large number of atoms, the time and computing capacity necessary is too large if we use the Fock basis. In fact, some 
data are missing for the Fock basis due the long computing times needed to generate them. The comparison between 
computing time is important in order to economize computing resources, because this permits to calculate systems 
with more atoms in less time, in other words for large N the coherent basis optimizes the use of the computing 
resources. It can be seen that the time and number of excitations N max needed for convergence decreases in a 
significant way if we use instead the coherent basis. 



In figure 3 and figure 4 we show the ground and first excited energies and in figure 5 the difference between them 
varying the interaction parameter. Likewise, in figures 6 to 9 it is shown the photon number and number of excitations 
expected values as functions of 7. Both observables are divided by j in order to scale them with a system natural 
size. It can be seen in all the observables that there is a critical value of the interaction parameter where a phase 
transition takes place. Also, we note the effect is clearer with the increase of j, but we do not need large values of 
N to observe the transition. So, using this basis, we can reproduce the results about the observables shown in [8[. 
As it was shown in @-[ll|, coherent states are a good asymptotic approximation for the states as j tends to infinity 
and could be used as trial functions to this kind of Hamiltonian. Therefore, it is not surprising that these coherent 
states let us to obtain easily a numerical description of the Dicke Hamiltonian. With these tool we intend to explore 
another observables of interest in the model. 
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For j=l 


"< 


Energy 


Fock Umax 


Fock Computing time 


Coherent N raax 


Coherent Computing time 





-1.00000 


2 


0.015 


2 


0.063 


0.1 


-1.00504 


3 


0.031 


3 


0.081 


0.2 


-1.02062 


4 


0.094 


4 


0.296 


0.5 


-1.15370 


7 


0.218 


8 


2.200 


1.0 


-2.15428 


13 


0.967 


16 


27.503 


2.0 


-8.03226 


27 


6.786 


37 


806.712 


For j=2 


"< 


Energy 


Fock n m ax 


Fock Computing time 


Coherent N raax 


Coherent Computing time 





-2.00000 


2 


0.011 


2 


0.01 


0.1 


-2.00505 


3 


0.093 


3 


0.265 


0.2 


-2.02085 


4 


0.171 


4 


0.499 


0.5 


-2.17336 


8 


0.702 


8 


4.430 


1.0 


-4.26487 


19 


6.692 


15 


46.426 


2.0 


-16.06350 


42 


65.520 


24 


296.495 


For j=5 


1 


Energy 


Fock Umax 


Fock Computing time 


Coherent N ma x 


Coherent Computing time 





-5.00000 


2 


0.218 


2 


0.702 


0.1 


-5.00506 


3 


0.437 


3 


1.888 


0.2 


-5.02099 


4 


0.717 


4 


3.196 


0.5 


-5.19716 


10 


7.052 


8 


35.007 


1.0 


-10.63840 


30 


150.354 


12 


144.395 


2.0 


-40.15720 


78 


2499.9600 


16 


425.243 


For j=10 





Energy 


Fock Umax 


Fock Computing time 


Coherent N mm 


Coherent Computing time 





-10.00000 


2 


0.686 


2 


2.137 


0.1 


-10.00510 


3 


1.529 


3 


5.554 


0.2 


-10.02100 


5 


4.446 


4 


11.481 


0.5 


-10.21300 


13 


50.919 


8 


102.836 


1.0 


-21.26310 


47 


2104.800 


10 


222.161 


2.0 


-80.31340 






13 


583.304 


For j=20 





Energy 


Fock Umax 


Fock Computing time 


Coherent Nmax 


Coherent Computing time 





-20.00000 


2 


2.824 


2 


8.206 


0.1 


-20.00510 


3 


5.647 


3 


22.105 


0.2 


-20.02110 


4 


17.831 


4 


46.785 


0.5 


-20.22690 


15 


311.003 


8 


424.572 


1.0 


-42.51290 






9 


635.408 


2.0 


-160.62600 






10 


888.784 



Table I: Computing times and n max (N max ) for several j as 7 is increased for both basis. 



IV. CONCLUSION 



Wc have compared the numerical solutions for the ground and the first excited states for the finite Dicke Hamiltonian 
employing the Fock basis and the coherent basis. In the case of the Fock basis it was shown that both the computing 
time and the number of excitations in the field necessary for the convergence of the solution increases as the number 
of atoms do. In contrast, for the coherent basis the photon number needed to obtain convergence diminishes as 
the interaction parameter grows. Also, the expectation values of the principal observables agree with the ones in 
Q. Finally, we have a tool which let us study more about the Hamiltonian, saving computing resources and tell us 
something about the system. 

JGH thanks B.M. Rodriguez-Lara for pointing out Ref. [3]. This work was partially supported by CONACyT- 
Mexico, FONCICYT (project-94142), and DGAPA-UNAM . 
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Figure 8: Expectation value for J z for the ground state as function of 7, for j from 1 to 20. 
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Figure 9: Expectation value for J z for the first excited state as function of 7, for j from 1 to 20. 
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